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DEVELOPMENT OF DIRECT-INVERSE 3-D METHODS 


FOR APPLIED TRANSONIC AERODYNAMIC DESIGN AND ANALYSIS 

I . Introduct ion 

This report cowers the period -from Janaury 1, 1987 thru June 30, 
1987. The primary tasks during this period were the completion of tasks 
associated with the -first phase of the project, the initiation of work 
involwing wiscous interaction effects, and the continued development and 
testing of design methods. 


1 1 . Personnel 

The staff assoc i ated wi th this project during the present reporting 
per i od were : 

LelandA. Carlson, Principal Investigator 
January thru May, Approximately 1/8 time 
June — Approximately 3/4 time 

Thomas Gaily, Graduate Research Assistant 
January thru June 

Robert Ratcliff, Graduate Research Assistant 
June — 

It should be noted that the first phase of the research work 
associated with this project has formed the basis for the Masters Thesis 
of Mr. Gaily, who received his M.S. degree in Aerospace Engineering in 
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May. It is planned to use the research in the second phase o-f the 
project as the basis -for the thesis o-f Mr. Ratcliff. 

III. Research Prooress 

As indicated above, the primary task during this reporting period 
has been to “complete" the -first phase o-f the project, which was to 
develop a versatile inviscid direct-inverse wing design method and 
program. The work in this area is presented in detail in Re-ference 1 
and is summarized in Re-ference 2. A copy o-f Ref. 1 has been sent to 
NASA and Ref. 2 is included as an appendix to this report. 

As part of this work, the methodology of handling surface relofting 
has been improved so that the method can now treat two difficult design 
tasks. The first task was to change a wing from super-critical to sub 
critical, which requires large changes in the leading edge region 

through relofting; and the second task was to make large surface changes 

to an orginal set of airfoil sections without generating errors due to a 
large number of geometry calculations. These objectives have been 

successfully demonstrated and are presented as "Test Case F" in the 

paper included in the Appendix. In this example, the wing for the 
Lockheed Wing-Body A planform at supercritical conditions was changed 
from a tapered thickness ranging from 127. to 67. to a 67 thick 
subcr-itical wing. Currently, efforts are in progress to go the other 
direction and design a supercritical thick wing starting from a thin 
subcritical set of airfoil sections. 
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Also, during this reporting period, work has begun on including 
viscous interaction effects in the design process. In the present 
context, it is planned to include as part of the viscous interaction 
e-f-fects the wing surface boundary layer, the wake thickness, and the 
wake curvature, either individually or together. Some preliminary 
viscous results are shown on Figures 1-4. 

Figure 1 shows the design situation, which is very similar to 
Design Case C discussed in the Appendix. For this case, two 

discontinuous wing sections were designed using pressures obtained from 
a viscous analysis of' a preselected wing planform. The wing sections 
consisted of NACA 0012 sections outside the design regions and modified 
NACA 0012 sections within the design regions. The pressures obtained 
from the analysis are shown in Figure 2 for a Reynolds number of 
approximately 11 million along with the inviscid analysis pressures used 
for Design Case C. The effect of including the viscous interaction can 
easily be seen in the more forward location and weaker shock strength 
for the viscous case, particularly at the outboard stations. 

Figure 3 compares the initial, target, and final design surfaces 
for this case. As can be seen, the target surfaces were accurately 
obtained for each section in the inverse region. However, slight 
deviations were present near the trailing edges on the stations 

bordering the direct-analysis regions. While this phenomena has been 
observed for inviscid cases and is believed to be due to spanwise slope 
variations, the deviation for this case was more pronounced. In 

addition, detailed examination of several inviscid cases has revealed 
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that sometimes- there appears to be an every-other spanwise point 
coupling in the designed air-foil sections. The origin of both of these 
phenomena is not de-finitiely known; and, consequently, both are 
currently being investigated. 

0-f course, the proo-f o-f a successful inverse design procedure rests 
in having an analysis of the designed wing predict the pressures 
originally desired. Such a comparison is shown on Figure 4; and, as can 
be seen, the agreement is quite good for preliminary design. 

Another effort which has been initiated during this reporting 
period has been the development of an additional design strategy option 
to the program. Currently, the method inputs pressures at arbitrary 
stations, linearly interpolates these to the appropriate computational 
grid points within the inverse region, and calculates inverse boundary 
conditions and surface changes at each of these points. The new option 
will specify inverse boundary conditions and will calculate surface 
changes (i.e. new airfoil sections) only along grid lines for which a 
pressure distribution has been input. If there are grid lines inbetween 
these two sections, then that region will be treated using analysis type 
of boundary conditions. However, each time the boundary airfoil 
sections are updated, the sections inbetween will be updated using 
linear interpolation. Thus, in this option, the user will only specify 
pressure distributions at the inboard and outboard stations of the 
design region; and the new airfoil sections in the design region will 
vary smoothly along with, hopefully, the resultant pressure 


distributions. 


While this, effort is still in the developmental stage and no 
results are yet available, i f it is success-ful it may offer several 
advantges. First, the new input ■format o-f inputting pressures at the 
computational grid lines, may ease the current sensitivity o-f the method 
to input pressure distributions, particularly when shock waves are 
included and when the present method linearly interpolates to obtain 
values at grid lines. Second, and perhaps most importantly, this option 
may be more applicable to many engineering situations and 
computationally ■faster. 

IV. Future Efforts 

During the next reporting period it is anticipated that most o-f 
the viscous interaction developmental work and new design option work 
will be completed. In addition, it is hoped that detailed verification 
studies can be started. 


V. Grant Mon i tor 

The NASA Technical Monitor for this project is Richard L. Campbell, 
Applied Aerodynamics Group, NTF Aerodynamics Branch, Transonic 
Aerodynamics Division, NASA Langley. 
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INVISC1D TRANSONIC WING DESIGN USING INVERSE 
METHODS IN CURVILINEAR COORDINATES 


Thomas A. Gaily* Leland A. Carlson** 

Texas A&M University Texas A&M University 

College Station, Texas College Station, Texas 


ab stract 

An inverse wing design method has been developed 
around an existing transonic wing analysis code. The 
original analysis code. TAW FIVE, has as its core the 
numerical potential flow solver, FLO SO, developed by 
Jameson and Caughey. Features of the analysis code 
include a finite -volume formulation; wing and fuselage 
fitted, curvilinear grid mesh ; and a viscous boundary 
layer correction that also accounts for viscous wake 
thickness and curvature. The development of the inverse 
methods as an extension of previous methods existing for 
design in Cartesian coordinates is presented. Results are 
shown for inviscid wing design cases in super-critical 
flow regimes. The test cases selected also demonstrate 
the versatility of the design method in designing an 
entire wing or discontinuous sections of a wing. 

m m^ QL^iuEL 

Cp - Coefficient of pressure 

h - Jacobian of coordinate transformation 

H - Jacobian matrix 

j - Transpose of inverse Jacobian matrix 

Mqo ~ Freestream Mach number 

q^ - Magnitude of freestream velocity 

Q - Magnitude of local velocity 

u,v,w - Components of physical velocity vector 

U,V,W - Components of contravariant velocity vector 

a - Angle of attack 

7 - Ratio of specific heats 

6 - Differential operator 

6(x) - Displacement thickness 

$ r (x) - Displacement thickness due to relofting 

A - Trailing edge thickness 

A t - User specified trailing edge thickness 

p - Averaging operator 

p - Density 

- Reduced/perturbation potential function 

($> - 4> + x cos(a) + y sin(a) ) 

& - Potential function 

INTRODUCTION 

In recent years the importance of transonic flight to 
both military and commercial aircraft and the develop- 
ment of specialized transonic wings for several flight 
research experiments have prompted significant efforts to 
develop accurate and reliable computational methods for 
the analysis and design of transonic wings. Many methods 
of solution have been developed, but among those which 
have shown promise due to their computational efficiency 
and engineering accuracy have been those based upon the 
full potential flow equations in either their conservative 
or non-conservative form**^. The TAWFIVE 4 FORTRAN 


• Graduate research assistant. 

** Professor of Aerospace Engineering, Associate 
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code in particular has proven to be an excellent and 
reliable analysis tool. This analysis code is based upon the 
FLO30 finite volume potential flow method that was 

developed by Jameson and Caughey^. Among the fea- 

tures of FLO30 are its fully conservative formulation and 
its three-dimensional curvilinear grid. The latter can be 
fit around any general combination of fuselage shape and 
wing planform. 

The purpose of the research described in this paper 
has been to develop a wing design method that is based 
on the existing TAWFIVE analysis code and is compatible 
with the existing computational methods and program 
structure of that code. Of the many wing and airfoil 
design methods available 5 * 8 , the inverse method as 
developed by Carlson^*^ was selected for use. The 

current work extends the previously developed design 

methods developed for orthogonal grids to the more 
generalized curvilinear grid system of TAWFIVE, while 
also providing greater design flexibility and versatility for 
engineering applications. These last goals were achieved 
by the inclusion of user options for designing either the 
entire wing or only discontinuous wing segments as 
shown in Figure 1. The availability of this option is 
useful to engineers who are typically faced with desig- 
ning around regions where the wing geometry may be 
fixed by constraints other than aerodynamic consider- 
ations. 
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WINC. ANALYSIS METHODS 
Potential Flow Solver 

The inviscid potential analysis of TAWFIVE is 
performed by the program FLO30 developed by Caughey 
and Jameson^’ For a complete description of the 
FLO30 code and its theoretical basis the reader is 
referred to Caughey and Jameson’s papers and some 
earlier developmental work by Jameson* 4 ~*^. A brief 
description is presented here to provide for completeness 
and to provide a background for the inverse design 
developments which will be discussed in detail. 

FLO30 solves the full potential equation in conserva- 
tive form which when transformed from Cartesian coor- 
dinates to generalized curvilinear coordinates is: 

(phU) e + (phV)jj + (phW) c = 0 (1) 

where the subscripts denote differentiation with respect to 
the curvilinear coordinates £, and f. The contravariant 
velocities are related to the physical velocities and the 
derivatives of the potential function by: 


0 - H “ " ,T “ rl If! <a 

and H is the transformation matrix defined by: 

[ V x„ x r l 

y £ y n y f wi,h h ” l H l ( 3 > 

z £ ^ z fJ 

The local density can be obtained from isentropic 
relations as: 

P » [1 + ^"5 Moo^ 1 - u 2 - v 2 - w 2 )] 7 1 (4) 

The numerical approach used in FLO30 is a finite 
volume technique. To understand this approach, consider 
the simple two dimensional case represented by the grid 
system shown in Figure 2. 
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Figure 2. Finite-Volume Cell Location 


The dashed cube shown in the figure indicates the area 
element under consideration. The flux of fluid through 
side a-b can be approximated by the average of the 
fluxes at point a and b with similar results for the side 
c-d. The net flux in the x direction for the elemental 
area centered at grid point i,j is then: 

(phU)£ - ((phU a + phU b ) - (phU c + phU d >] / 2A£ 
or in the notation of Caughey and Jameson, 

(phu)^ * *y € <pu) 


where /i indicates averaging and 6 indicates differentiation 
in the indicated directions which are defined as follows 
(allowing A£=Ar?=Af=I): 


' C '* ,J 1 — j , j , K 

(^U)i ; k = ( u i+i,j,k + u i-i,j,k >/ 2 

,i,k “ ^i+i,j+i,k + , j -i ,k + ^M,j+l.k 


. . etc. 


When extended to the other flux components and to 
averaging over cube surfaces in three dimensions, the 
numerical potential equation is of the form: 

t*r}<; s £(P hV ) + ^y/> hv ) + py f (phW) * 0 

To find the flux quantities phU, phV, and phW at 
the finite volume cell vertices (i.e. points a, b, c, and d 
for the two dimensional case), it is necessary to evaluate 
Equations (2) through (4). The derivatives in these 
expressions can be expanded by the same volume averag- 
ing approach used above, thus: 


In m ¥( s v? 


x £ = 

y c = ^c¥ y > 
z £ “ 


with similar terms for the other transformation metrics. 
The above expressions, being centered at grid midpoints, 
will involve the values of the potential and grid position 
at grid points which are known from the previous poten- 
tial solution and the grid geometry, respectively. 


When solving transonic flows it is necessary to 
include in the solution algorithm some form of supersonic 
upstream dependence in order to account for both the 
physical nature of the flow and the viscous nature of 
shock waves, respectively. Caughey and Jameson intro- 
duced upwinding by the addition of terms into their 
potential numerical equation which are only non-zero 
when the flow is supersonic. Also, the finite volume 
technique exhibits a tendency for uncoupling of the flow 
field solution between alternating grid points. As a 

result, additional terms are included in the numerical 
potential equation. The final numerical equation which is 
solved by FLO30 when these terms have been included 
has the form: 




where P, Q, and R are the upwinding terms and Q^, 
Q?#’ anc * Qfrtf are the decoupling terms. 


Computational Grid Geometry 


The computational grid used by FLO30 is a body 
fitted, non-orthogonal curvilinear mesh constructed about 
a wing fuselage combination. The number of grid points 
composing the computational domain is typically 40 x 6 x 
8, 80 x 12 x 16, or 160 x 24 x 32 for the number of 
77, and f points in the coarse, medium, and fine grids, 
respectively. The grid is conformally mapped to the 
wing and fuselage surfaces as can be seen from the plot 
of surface grid lines shown in Figure 3. 

The grid is formed around spanwise airfoil sections 
in a similar manner in which "C" grids are mapped to 
airfoils in two-dimensional analysis. In addition, each 
spanwise computational plane is also conformally wrapped 
about the fuselage surface and a line extending forward 
from the fuselage nose. 
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Figure 3. Surface Grid Point Geometry 

A final set of grid surfaces are generated beneath 
the wing and fuselage surfaces and beyond the symmetric 
plane in order to aid in the formulation of both the 
finite-volume numerical flow equations and the flow 
tangency boundary conditions upon these boundaries. 
The grid points composing the "ghost" surfaces are 
calculated from linear extrapolations of the computation 
grid lines from inside the physical domain. 


conditions which need to be imposed at the slit is that 
the flow velocities, and thus pressure, be continuous 
across the cut. The flow potential, however, will have a 
discontinuous jump across the wake which is proportional 
to the sectional wing lift coefficient. 

INVERSE WING DESIGN METHODS 

As stated previously, a direct-inverse approach to 
wing design was selected for incorporation into the 
TAWFIVE code. The direct-inverse method derives its 
name from the division of the design wing surface into a 
fixed geometry leading edge region, where flow tangency 
boundary conditions are imposed, and an aft, variable 
geometry section where pressure boundary conditions are 
enforced. The pressure boundary where the user speci- 
fied pressure distributions are imposed does not extend 
forward to the leading edge due to difficulties of 
enforcing this type boundary condition near the beginning 
of an airfoil section. This restriction on the size of the 
pressure specification region does not seriously reduce the 
versatility of the design method since the leading edge 
regions for most airfoils are similar, and it is relatively 
easy to select a leading edge geometry which will 
produce the desired Mach number or pressure values at 
the beginning of the inverse region. In addition, specific 
leading edge shapes may be required due to other design 
constraints such as the necessity to house a leading edge 
flap or slot system. 

Press ur e , Bp u n .(jar.y.,£g,ndttifln 


Boundary Conditions 

Since the governing potential equations are written in 
terms of perturbations from free-stream conditions, the 
subsonic, far-field requirement that the flow return to 
the free-stream velocity and direction is satisfied by 
setting the perturbation potential equal to zero on the 
side and upstream boundaries. The downstream boundary 
condition is a "zero" order extrapolation of the potential 
(constant potential assumption) to the outflow boundaries. 

A flow tangency condition is applied along both the 
wing and fuselage solid surfaces by setting the normal 
contravariant component of the velocity vector to zero on 
the surfaces. This condition provides an equation which 
when approximated by a finite-difference expansion 
about the surface grid points can be used to set a value 
for the perturbation potential on the "ghost" grid points 
beiow each surface. Note that this finite-difference 
boundary condition differs in formulation from the 
finite-volume solution algorithm of the governing 
equations. As a result, it is possible to impose flow 
tangency using the finite-difference technique yet still 
have a slight normal surface velocity when performing 
the finite-volume calculations. Since it is essential to 
have accurate' boundary conditions at the wing surface in 
order to generate accurate solutions, a second condition is 
imposed upon the wing surface. This additional condi- 
tion involves reflecting the flux quantities calculated by 
the flow solver for the cell centers directly above the 
wing surface to the "ghost" cell centers beneath. The 
reflected normal fluxes then cancel each other out in the 
residual expression and a net zero flow is obtained 
through the surface. Similarly, a zero flux condition is 
applied at the half body symmetric plane, limiting 
solutions to symmetric, non-sideslip cases. 

The trailing edge slit boundary is not an actual limit 
to the physical domain as the other boundaries are, but is 
simply an artificial boundary created by unwrapping the 
physical plane into the Computational domain. The only 


In the inverse design regions on the wing, a pressure 
boundary condition will be specified rather than the flow 
tangency condition used in analysis zones. In formulating 
this boundary condition it is necessary to relate the user 
specified pressure coefficient, C p , to the current 
perturbation potentials at inverse design grid points. 
Consider the full potential equation for the pressure 
coefficient: 


C 


P 


- — 2 { [> M ood 



where: Q 2 


7 7 7 


If it is assumed that the pressure coefficient is 
primarily a function of the chord wise component of the 
velocity, u, and only slightly affected by the vertical and 
span wise components of velocity, v and w, then a stable 
approximation is made by time lagging the latter two 
velocities in the boundary condition expression. This 
assumption is true everywhere except near the leading 
edge; but since the inverse design boundaries have 
already been restricted to regions well behind the lending 
edge, the simplification is justified. The value of the 
local velocity, u, can then be calculated from the above 
expression in terms of the desired pressure coefficient 
and the current values for the vertical and spanwise 
velocities. In addition, the velocity u can also be 
calculated from the perturbation potentials using the 
relations of Eq. (2). Defining Jjj to be the elements of 
the inverse transpose of the Jacobian matrix, H, the two 
equations for u yield: 


Jll^+J|2^ + J|3^f “ Th 1 


1 2 2 [/ 

2 1 T 

' iMooCpX 


} * 2 )- >J 


i ♦ (S) 2 ♦ (g) 2 


- cos(a) 


(5) 


3 


ORIGINAL PAGE 13 

OF POOR QUALITY 


Since the spanwise and vertical flow velocities have 
already been assumed to be constant in the boundary 
condition relation, it is consistent to make the same 
approximation in the above expression with respect to the 
spanwise and vertical derivative terms, ^ and <t>^. This 
assumption is similar to the previous one, and leads to an 
explicit expression for the potential at one point. 

The finite difference approximation used involves 
expanding the derivatives of the potential about the 
mid-point i-i,j,k. The £ derivative is determined by a 
central difference involving the preceding and following 
grid point values. The rj and f derivatives are found at 
the mid- point by averaging the derivatives from the 
preceding and following grid points found by a three 
point backwards and central difference approximations, 
respectively. Figure 4 shows the point dependence and 
pressure specification point for this method. The 
resulting numerical expression obtained with these finite 
approximations is: 


n+J n 

r n+1 n n n 

♦ J l2L?<*i,j,k + *i-l,j,k) “ 4 (*i,j-l,k + *i-l j-I t k> 
n n -» 

+ *i,j-2,k + ^i-l,j-2,kj t 4 
n n n n 

+ J 13<*i.j,k+1 + *i-I,j,k+l " *ij,k-l * <*i-I,j,k-I >/ 4 

- F(CPi_* >k ) 


Here, the superscripts n and n+1 refer to current 
values of the potential and the new values of the 
potential being imposed by the boundary condition, 
respectively. Also, the term F(Cpi_$ |<) is the right hand 
side of Eq. (5) evaluated using the' pressure coefficient 
specified at point M,k. Solving the above expression for 
the potential at point i,j,k yields: 


n+1 

^i.j.k 


J n+ 3J 12 /4 


1 J ll*i-l,j,k 

1 n n n 

- J^L^i-lJ.k * 4 ^i,j-l,k + +i-I.j-l.k> 
n n , l 

+ ^i,j-2,k + ^i-l,j-2,kj / 4 
n n n n 

J 13Wi,j,k+l + ^i-l,j,k+l * ^i,j,k-l " 1 ,j,k- 1 >/ 4 


+ p ( c Pi-i,k) 


The potential values at n+! in the direct region are 
known initially since they do not change when the 
inverse boundary condition is applied; i.e. ^ n+1 « 4> n . All 
the potentials on the inverse boundary can then be 
calculat d and, since the spanwise and vertical derivatives 
are small, will primarily be functions of the pressure 
coefficient at grid point i-* and the value of the 
potential at grid point i-1. 



Figure 4. Point Dependance and Location 



Surface Calculations 

As the inverse boundary conditions drive the flow 
field to a converged solution, it is necessary to 
periodically calculate the location of the new displacement 
surface and to regenerate the computational grid about 
this new geometry so that the pressure boundary surface 
will correspond to the physical boundary surface. Each 
new surface can be found relative to the previous surface 
from an integration of the wing surface slopes. However, 
the surface slopes must first be calculated from the 
current flow field solution using the flow tangency 
boundary condition which in curvilinear coordinates is: 

V T X VF = 0 


The only concern with using this mid-point specifi- 
cation scheme is that the current method of calculating 
the pressure data output from FLO30 uses a grid point 
centered difference scheme for the streamwise derivative. 
This difference could potentially allow a pressure to be 
specified correctly but still have a significantly different 
value output from FLO30 due to the inconsistent calcula- 
tion methods. However, as shown on Figure 5, where the 
pressures calculated for a typical flow solution are 
compared for the two different calculation techniques, 
this possible error has not been significant in practice. 


where V is the contravariant velocity vector and VF is 
the gradient of the surface function with respect to the 
curvilinear coordinates. Note this condition is a direct 
analog to the same condition expressed in physical space. 

A more useful expression can be obtained by 
expanding the above equation to: 



This expression can be solved for the new chordwise 
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airfoil slopes, drj/dZ, if the current values of the 
spanwise slope, dr\/d$ % are used. Since the wing surface 
is represented in the computational grid as a plane of 
constant rj, the current slopes on the wing surface equal 
zero and a simplified flow tangency condition results: 



The above expression has been applied to the com- 
putational surface plane in order to find the relative 
location of the new physical surface. This approach is 
an approximation since the above equation is only exactly 
true when applied to the new surface itself. Using this 
method, however, provides for a stable iterative surface 
updating procedure which quickly converge to the target 
surface. 

To calculate the relative surface slopes, it is first 
necessary to accurately determine the values of the 
contravariant velocities, U and V. As was also deter- 
mined by the work of Weed, et al. 12 , a simple finite 
difference calculation of these velocities is not 

sufficiently accurate. Borrowing from Weed, et a!., a 
more accurate method was implemented which uses the 
residual expression to calculate the velocity ratio, V/U, 
under the assumption that the residual is zero at the 
surface points. The residual expression from FLO30 can 
be written in finite volume form as: 

+ ^5 f (phW) 

♦ (other terms) * 0 

The "other terms* in the above expression involve the 
grid point coupling and upwind dependence terms of the 
formulation and are assumed to be constants in the 
following development. 

The desired velocities can also be written in this 
finite volume form as: 

V * phV - p^(phV) and U « phU « n^iphU) 

By simple manipulations, the normal velocity can be 
obtained from the residual expression as: 

2/^ c (phV) « 2^(phV) r/ ., - v n MphU) (6) 

- HfySfiphvf) - (other terms) 

where the subscript tj - ! refers to the values at grid cell 
centers above ihe wing surface. 

In order to use Eq. (6) to find the desired surface 
velocity ratio, it is necessary to know the U and W 
yelocity components at the "ghost* cell centers below the 
wing surface. These values can be obtained in a manner 
consistent with FLO30 by specifying the "ghost" cell 
values to equal the values at corresponding points 
immediately above the wring surface. A comparison of 
the accuracy of both the finite difference approach and 
residual approach is shown in Figure 6. The calculated 
displacements are for a converged analysis solution for 
which the calculated slopes should of course be zero. 

With the contravariant velocities known, an integra- 
tion of Eq. (6) through the inverse design region from 
the leading edge to the trailing edge yields a set of 
surface displacements, 6(x) % for the new wing surface 
relative to the previous one. These displacements are 
expressed as changes in the computational coordinate rj, 
and are converted to surface displacements in the 
physical plane via the, local grid transformation. The 
physical plane displacements are coincident with the 



X/C 


Figure 6. Comparison of Slope Calculation Methods 

computational grid points in the inverse regions. To 
obtain the corresponding displacements at the original 
geometrical locations specified in the program input data, 
a linear interpolation of the above data is performed. 
Finding the displacements at the original geometry 
stations permits the calculation of the new wing airfoil 
sections at the same semispan locations. 

Trailing Edge Closure 

The procedures outlined above will compute a wing 
surface corresponding to a given, fixed, leading edge 
geometry and to a desired set of pressure distributions in 
the inverse regions. The above procedures do not, 

however, guarantee that this wing geometry will be 
practical. In particular, past experience 9 has shown that 
inverse surface calculations may yield airfoil sections 
which have either excessively blunt trailing edges or 
which, at least numerically, have the upper and lower 
surfaces crossed at the trailing edge ("fish tailed”). The 
former case is undesirable due to aerodynamic consider- 
ations, while the latter is physically impossible and may 
produce unpredictable problems in the grid generation or 
flow calculation portions of FLO30. 

Since for any specified pressure distribution the 
corresponding wing surface will be controlled by the 
leading edge geometry, which serves as an initial spatial 
boundary condition for the inverse region, the problem of 
assuring trailing edge closure can be viewed as the proper 
seieciion of a leading edge shape. A procedure for 
systematically modifying the leading edge region in order 
to achieve some desired trailing edge thickness is called 
relofting. Such a relofting procedure has been incorpor- 
ated into the present design process in order to both 
prevent the problems of trailing edge crossover and to 
allow the user the option of specifying a trailing edge 
thickness as an additional design variable. This design 
feature should be very useful in practical applications 
since it automates the iterative selection of a leading edge 
shape which would otherwise have to be performed by 
the user. 

Two methods of relofting can presently be selected. 
The first method is a simple linear rotation scheme. This 
method can be visualized with the help of Figure 7. The 
dashed line indicates the original leading edge geometry 
and a hypothetical new surface shape which has been 
calculated for the inverse design regions. Without 
modification, this new surface has a trailing edge 
thickness of A. If a thickness of A t were specified by 
the user, then the surface would have to be relofted or 
changed. In the present scheme, in order to obtain the 
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Figure 7. Relofting to Force Trailing Edge Closure 

desired thickness, a displacement thickness, 5 r , is added 
to the current design surface. This thickness has a 
distribution from the leading to the trailing edge and is 
determined by the formula: 

S r (x) « (A t - A) (x/c) 

where c is the chord length of the local airfoil section. 
The total displacement for a surface update is then the 
sum of the two displacements, £(x) and £ r (x). When 
both the upper and lower surfaces are designed simulta- 
neously, the displacement magnitudes determined by 
relofting are divided between the two surfaces so that 
half is added to the lower surface and half to the upper 
surface. 

The second relofting method uses the same approach 
as the first for the aft inverse regions, but modifies the 
leading edge region by a proportional thining or thicken- 
ing of the surface ordinates. This approach can be 
expressed by: 


y n+1 U) 


V .n+I 
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where the j subscript refers to the ordinate at the 
direct-inverse junction determined from the linear 
relofting of the aft regions. Note that this method will 
produce leading edges in the same family of shapes and, 
for example, allow the design from a NACA 0012 airfoil 
to a NACA 0006 airfoil (see Test Case F). 


RESULTS 


A variety of different test cases were run as 
verification of the current design method. These cases 
involved both subcritical design and supercritical design 
over section geometries selected to test the versatility of 
the input and design control logic. In this section results 
from three of the more significant test cases will be 
presented. The results shown were obtained on a 
medium grid having 81 streamwise, 13 vertical, and 19 
spanwise points with 1 1 spanwise stations and 53 points 
on the wing at each station; and in all cases the 
maximum change in the reduced potential was reduced at 
least three orders of magnitude. Thus, the results do not 
represent ultimate convergence but should be represent- 
ative of '’engineering accuracy". 

The planform selected for the test cases was the 
Lockheed Wing A wing-body. The wing for this config- 
uration has a quarter chord sweep of 25 deg., a linear 
twist distribution ranging from 2.28 deg. at the wing 
body junction to -2.04 deg. at the wing tip, an aspect 
ratio of eight, and a taper ratio of 0.4. The last two 
values are based upon the wing without fuselage. 
However, instead of the supercritical sections normally 
associated with Wing A, the initial airfoil sections at each 
span station were assumed to be composed of symmetric 
NACA four digit airfoil sections. 

The target pressure distributions used in the design 


regions for the first two test cases were selected to yield 
airfoil shapes thicker in the aft portions of each section; 
and, at supercritical conditions, to yield on the upper 
surface weaker and more forward shock waves than those 
which would normally occur on a NACA 0012 section. 
On the lower surface, the target pressure distributions 
were selected to have either a favorable pressure gradient 
or fairly constant pressure plateau over much of the 
lower surface. 

For the last test case, the pressure distribution was 
obtained from analysis solutions of an assumed wing 
geometry. The intent of this cases is to verify the 
relofting procedures and show the ability of the current 
method to make large surface changes in going from a 
thick wing to a thin wing (approximately 12 percent to 6 
percent thick respectively). 

All cases were for a freestream Mach number of 0.8 
and an angle of attack of two degrees. In each case, the 
pressure distribution was specified in the design regions 
from the 15% local chord location to the trailing edge 
and used as the boundary condition in these inverse 
regions starting with the first iteration. Normally, three 
hundred SLOR iterations were executed prior to the first 
design surface update calculation; and subsequently, 
surface updates were computed every fifty cycles. 
Usually, the solution was considered converged and 
terminated after 550 total iterations for the first two 
cases and, due to the large amount of relofting required, 
after 950 iterations for the last case. 

T g st .Ca ss . Q 

The inverse design regions for Case C, which was an 
attempt to design both upper and lower surfaces on two 
noncontiguous regions of the wing at supercritical 
conditions, are shown on Figure 8; and a comparison 
between the initial pressure distribution associated with 
NACA 0012 sections and the target pressures for two 
sections is portrayed on Figure 9. As can be seen, the 
target pressure distribution essentially eliminates at 
inboard stations the upper surface shock wave present on 
the original wing; and at outboard stations it weakens the 
shock and moves it forward. In addition, significant 
changes in the lower surface pressure gradients are 
evident. Also shown on Figure 9 are the pressures 
computed by the program at the end of the inverse 
design procedure (denoted as "design pressures"). These 
pressures are in excellent agreement with the target 
pressures, which indicates that the method is satisfying 
properly the desired inverse boundary conditions. 

The corresponding designed airfoil sections for this 
case are shown on Figure 10. Even on the expanded 
scale, the agreement between the designed and target 
surfaces is excellent at all design stations. However, 
trailing edge closure was not enforced for this case; and 
there is at the boundary stations some departure between 
the designed surfaces and the target surfaces near the 
trailing edge. It is believed that this slight difference is 
a ramification of the change in spanwise slopes near the 
trailing edge between the direct and inverse regions. 

In any event, the pressure distributions resulting 
from an analysis of the designed surfaces shown in 
Figure 10 are in excellent agreement with the target 
pressures, as can be seen on Figure II. In addition, the 
section lift coefficients at the various design stations are 
in very good agreement with the target coeft icients. 
Based upon these results it is believed that the present 
method can adequately design/modify nonadjacent regions 
of a wing in transonic flow. 
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Figure 8. Design Case C 



Figure 9. Comparison of Initial Pressures with 
Target Values (Case C) 



Figure 10. Comparison of Designed Sections with 
Original and Target Sections (Case C) 




Figure II. Comparison of Pressures from Analysis of 
Designed Wing Pressure Distributions (Case C) 


Test Case E 

For this test case, it was decided to design two 
non-adjacent upper surface regions simultaneously with a 
lower surface region which overlapped the upper zones. 
The location of these inverse design regions is shown on 
Figure 12. Likewise, Figure 13 compares the pressures 

associated with the initial wing sections shapes to the 
target pressures and to the pressures computed at the end 
of the design calculation for three design stations. It 
should be noted that this case is for supercritical 
condition and trailing edge closure is not enforced. As 
can be seen, at stations where only one surface is being 
designed (e.g. 50%, and 70%) the pressure distribution on 
the fixed surface also changes due to three dimensional 
effects from adjacent station which have been redesigned. 
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However, as depicted on Figure 14, only the design 
surfaces change form the original shape; and these 
surfaces are in reasonable agreement with the target 
profiles. 

Finally, Figure 15 compares analysis results obtained 
for the designed wing with the target pressures. Even 
for this complicated case, the agreement between the two 
distributions and between the actual and target lift 
coefficients is excellent. 



Figure 12. Design Case E 



Figure 13. Comparison of Initial Pressures with 
Targei Values (Case E) 



Figure 13. Continued 
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Figure 15. Comparison of Pressures from Analysis of 
Designed Wing Target Distributions (Case E) 


TFST CASE F 

The final test case was selected to demonstrate the 
ability of the design methodology to handle two difficult 
design tasks. The first task was to change a wing from 
super-critical to sub-critical. Due to the upstream 
dependance of the supersonic flow, this required making 
large changes in the leading edge region through the 
relofting procedures. The second task was to make large 
surface changes to the original airfoil without generating 
large surface distortions from the accumulation of 
geometry calculation errors. The design regions for this 
case are shown in Figure 16 where the wing thickness 
varied from 12% to 6% between the wing root and 80% 
span location and was constant going outward to the tip. 
The input design pressivres were for a constant 6% thick 
wing. 


The first attempts at this design used the linear 
leading edge relofting procedure and from a practical 
standpoint were unsuccessful. The final design surfaces 
were still supersonic in the leading edge regions while 
satisfying the subsonic aft surface conditions by 
producing strong shocks at the direct-inverse junction 
location. In addition, the surfaces themselves had sharp 
surface slope discontinuities at the same location. 

When the thining approach was used to reloft the 
leading edge, much better solutions were obtained. 
Figures 17 through 19 show the changes in pressure 
distribution and surface shapes with a comparison of 
target to designed surface pressures for a few span 
sections as in the previous cases. As can be seen, 
excellent agreement between target and final pressures 
and surface were again attained for this extreme case. 
The only noticeable surface irregularities are a small 
wiggle at the direct-inverse junction which can also be 
seen as a small pressure jump in Figures 17 and 19. 



Figure 16. Design Case F 



Figure 17. Comparison of Initial Pressures with 
Target Values (Case F) 
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Figure 17. Continued 


Figure 18. Continued 



Figure IS. Comparison of Designed Sections with 
Original and Target Sections (Case F) 



Figure 19. Comparison of Pressures from Analysis of 
Designed Wing Target Distributions (Case F) 
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CONCLUSIONS AND SUGGESTIONS 
EQR FUTURE WORK 

A direct- inverse wing design method has been suc- 
cessfully incorporated into the TAWFIVE transonic 
wing- body analysis computer code. The resultant code is 
capable of designing or modifying wings at both tran- 
sonic and subsonic conditions and includes the effects of 
wing- body interactions. A series of test cases have been 
presented which demonstrate the accuracy and versatility 
of this inverse method. 

Inclusion of viscous effects via the addition of the 
wing surface displacement thickness and wake thickness 
when performing wing design has been accomplished but 
not completely verified. Additional work will be 
required to run a sufficient sampling of test cases for 
evaluation of this design mode. The unique problems 
associated with viscous design and the effects of the 
various viscous correction models available in TAWFIVE 
would be the subject of a continuing research effort. 

The development and evaluation of alternate methods 
of surface relofting are also topics for which continued 
research is suggested. The current method of relofting 
restricts the user to a family of leading edge geometries 
which can be constructed by the linear rotation of the 
initial shape. The option of using other relofting 
methods would extend the family of available shapes and 
add versatility to the design method. 
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